Numerical inversion of a general incomplete elliptic integral

We present a numerical method to invert a general incomplete elliptic integral with respect to its argument and/or amplitude. The method obtains a solution by bisection accelerated by the half argument formulas and the addition theorems to evaluate the incomplete elliptic integrals and Jacobian elliptic functions required in the course. If a faster execution is desirable at the cost of complexity of algorithm, the sequence of bisection is switched on the way to the improvement by the Newton method, Halley’s method, or higher order Schröder methods. In the improvement process, the elliptic integrals and functions are computed by using Maclaurin series expansion and addition theorems based on the values obtained at the end of bisection. Also the derivatives of the elliptic integrals and functions are recursively evaluated from their values. By adopting 0.2 as the critical value of the length of solution interval to shift to the improvement process, we suppress the expected number of bisections as low as 4 in average. The typical number of application of update formulas in the double precision environment is three for the Newton method, and two for Halley’s method or higher order Schröder methods. Whether the improvement process is added or not, our method requires none of procedures to compute the incomplete elliptic integrals and Jacobian elliptic functions but those to evaluate the complete elliptic integrals once at the beginning. As a result, it runs fairly fast in general. For example, when using the improvement process, it is around 2-5 times faster than the Newton method using Boyd’s starter [3] in inverting E(φ|m), Legendre’s incomplete elliptic integral of the second kind. Email address: [email protected] (Toshio Fukushima) Preprint submitted to Journal of Computational and Applied MathematicsNovember 27, 2012